Universal behavior of relaxational heterogeneity in glasses and liquids 
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We report an investigation of the heterogeneity in super-cooled liquids and glasses using the non- 
Gaussianity parameter. We simulate selenium and a binary Lennard- Jones system by molecular 
dynamics. In the non-Gaussianity three time domains can be distinguished. First there is an 
increase on the ps-scale due to the vibrational (ballistic) motion of the atoms. This is followed, 
on an intermediate time-scale, by a growth, due to local relaxations (/3-relaxation) at not too high 
temperatures. A maximum is reached at times corresponding to long range diffusion (a-relaxation). 
At long times the non-Gaussianity slowly drops, the system becoming homogeneous on these time- 
scales. In both systems studied, the non-Gaussianity follows in the intermediate time domain, 
corresponding to the /3-relaxations, a law oc This general behavior is explained by collective 
hopping and dynamic heterogeneity. We support this finding by a model calculation. 



PACS numbers: 61.20.Lc, 61.43.Fs, 64.70.Pf 



Although glass is one of the most common materi- 
als, its physics and especially its dynamics are still only 
poorly understood. In addition to sound waves, two level 
systems 0J^] and quasi-local (resonant) vibrations jj), 
experiments indicate a wide distribution of relaxations, 
i.e. non-periodic changes of the local structure In 
super-cooled liquids one observes, apart from the vibra- 
tional (ballistic) motion of the atoms, two types of re- 
laxations with different time-scales. These [3- and a- 
relaxations are attributed to short range (cage) motion 
and diffusion, respectively. 

In recent years these relaxations, both in glasses and 
in liquids, have been studied intensively by experiment 
|^-^| and theoretically ]l0| |l3|| . One particular aim was 
to determine whether the relaxations involve only groups 
of atoms or are spread over the whole system. The first 
case, where relaxations are restricted to a few atoms only, 
is known as heterogeneous scenario; the other one as ho- 
mogeneous. 

Spatial heterogeneity is thought to be responsible for 
the non-exponential relaxations in super-cooled liquids 
. This view has recently been challenged on the basis 
of inelastic neutron scattering experiments on polymers 



however, see also 15 



To understand its effects, it is necessary to know the 
properties of the "dynamic heterogeneity" itself, e.g. the 
time and temperature dependence. Qualitatively it is 
known |i 12 1(|] , that the system becomes homogeneous 
at all temperatures for sufficiently long times, corre- 
sponding to the a-regime. In the intermediate time do- 
main, corresponding to the /3-relaxation, heterogeneity 
becomes more pronounced when the system is cooled 
down. Here we will show that there is a universal law 
governing heterogeneity at these time-scales. 

This is closely related to collectivity of motion. Mea- 
surements of the isotope effect have shown that diffusion 
both in glasses and in super-cooled liquids is highly col- 
lective |l7],[l8). A similar very small isotope effect was 
also observed in simulations of a Lennard- Jones liquid 



p9[ . In glasses, one observes collective jumps of chain- 
like structures [pO| [2^1 . Similar mobile structures are also 
observed in the under-cooled liquid p3"|]24| ]. 

Following our previous work |l6|] , we investigate the 
non-Gaussianity behavior, i.e. the heterogeneity, of re- 
laxations in glasses. In this letter we focus especially 
on intermediate times, shorter than the typical diffusion 
time. As shown previously, the non-Gaussianity increases 
markedly in this time domain. Here we want to go one 
step further. First we show that molecular dynamics sim- 
ulations of two different systems, Se and binary Lennard- 
Jones (LJ), give strikingly the same law for the increase 
of the non-Gaussianity parameter (NGP), av{t) defined 
further down. In this intermediate time domain the non- 
Gaussianity follows a power law a%(t) oc y/i, for both 
systems and both for temperatures above and below T g . 
We propose a simple model based on the collective hop- 
ping of groups of particles. 

The simulations for Se and LJ were both done with a 
velocity- Verlet algorithm, controlling the temperature by 
velocity adjustment and using the equilibrium volume at 
the given temperature, i. e. zero average pressure. 

We describe Se with a 3-body potential This po- 
tential has been used previously to calculate vibrations 
p6| , local relaxations [ p2|j23| and heat transport p7|] in 
amorphous Se. It provides a sound basis for the study 
of both the atomistic and the electronic structure p8[ . 
We prepared 4 independent samples of hot liquid, each 
containing 2000 atoms. These were then quenched to the 
desired temperatures with rates of 10 13 K/s. Before us- 
ing the configurations for the measurements, they were 
aged for several ns. The effective quench rates were thus 
of order 10 10 K/s. The glass transition temperature is 
estimated as T g w 300 K, and the critical temperature 
T c fa 330 K. More details are given in Ref. @. 

For the binary LJ simulations we take the frequently 
used parameters of Kob and Anderson [^2|,^9) . The sim- 
ulations are done with 5488 atoms and a composition of 
20% small particles. We quench from appropriately aged 
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samples with a rate of about 10 11 K/s (relating the LJ- 
values to Ar) and subsequently age the samples again. A 
heat bath was simulated by comparing the temperature 
averaged over 20 time steps with the nominal temper- 
ature. At each step 1% of the temperature difference 
was adjusted by random additions to the particle veloc- 
ities. To ease comparison with previous LJ simulations 
we shall, in the following, use the usual LJ units. De- 
tails of the simulation procedure are analogous to the 
ones described for the monatomic LJ system From 
the diffusion constant we determine T c k, 0.37. The dis- 
crepancy of this value compared to the one in previous 
simulations p| , p9[ , is due to the difference of densities. 
The previous work was done for a constant density of 
p = 1.2 whereas we find for our equilibrium samples at 
T = T c a lower value of p = 1.15. We also used a slightly 
larger cut-off of the potentials. 
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FIG. 1. Non-Gaussianity parameter multiplied by time 
against time in a log-log representation. The values are ob- 
tained from a molecular dynamics simulation of Se at the 
temperatures (from bottom to top): 495 K, 445 K, 400 K, 
355 K, 330 K, 290 K, 255 K and 200 K. 

To quantify the hetero geneity o f the relaxations we 

and use the non- 
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follow previous work 
Gaussianity parameter (NGP) 




a 2 (t) = - 



3 < Ar 4 (i) > 



5 < Ar 2 (t) > ; 



where < ... > denotes time averaging, Ar 2 (t) is the mean 
square displacement and Ar 4 (t) is the mean quartic dis- 
placement. Experimentally the NGP can be determined 
from the g-dependence of the Debye- Waller factor 0. 
From simulations the qualitative behavior of a 2 (t) is well 
known. Starting from a 2 {t = 0) = 0, it rises on a time- 
scale typical for vibrations (i w 1 ps) to values around 
a 2 = 0.2. In a hot liquid the NGP drops from this 
value on a ps-scale. In under-cooled liquids and in glasses 
the NGP keeps growing on intermediate time-scales and 
reaches values an order of magnitude larger. Only on 
the time-scale of diffusion or a-relaxation does the NGP 
drop and finally reaches, for t — > oo, the limit a 2 = 0. 



This latter limit reflects the ergodicity of the system for 
long times. From the increase, one concludes that the 
relaxations are mainly heterogeneous in the intermediate 
time scale. This becomes more and more pronounced as 
the system is cooled down [M . 

Here we want to stress a material independent prop- 
erty in this time domain. For this, we plot for Se a 2 (t) 
multiplied by t against time in a log-log representation, 
Fig. 1. The most interesting feature of this plot is the ap- 
pearance of a master curve for all temperatures stretching 
over a time domain from 10 _1 ps to 10 3 ps, i.e. four or- 
ders of magnitude. This master curve is independent of 
any rescaling constants or procedures. It corresponds to 
a power law i 3 / 2 leading to a 2 (t) oc \/t. 

To check whether this behavior results from the par- 
ticular structure of amorphous Se, which is constituted 
from chains and rings, we repeated the calculations for a 
binary LJ system. This model is frequently used as an 
idealized dense packed metallic glass. The nearest neigh- 
bor coordination is near 12 rather than 2 in Se. We plot 
the NGP of binary LJ in the same way as for Se, i.e. 
t ■ a 2 {t) versus t in a log-log plot. As Fig. 2 shows, the 
time dependence of the NGP follows the same power law 
a 2 (t) oc \fi as seen in Se. Moreover, we observe this be- 
havior not only for the average NGP, shown in Fig. 2, 
but also for both components separately. 




10 10 

t ( units of (m cr / e) ) 



10 



FIG. 2. Non-Gaussianity parameter multiplied by time ver- 
sus time obtained from molecular dynamics simulation of a 
binary Lennard- Jones system at the temperatures (from bot- 
tom to top): 0.88, 0.56, 0.48, 0.40, 0.36 and 0.32. 

From the above we conclude that the NGPs of differ- 
ent types of structural glasses and super-cooled liquids 
follow at intermediate times the same time dependence: 
a 2 (t) oc \fi. Therefore, we think that the mechanism 
responsible for the increase of the NGP, i.e. of the het- 
erogeneity, is common to all kinds of glass-forming ma- 
terials. 

In our previous investigation of the non-Gaussianity 
]l6[ , we have clearly shown that the increase of non- 
Gaussianity is due to relaxations. Moreover it has also 
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been shown that in under-cooled liquids and in the gl assy 
phase, clusters of so called mobile particles exist [[2l|-f24[ . 
Theses move in a given time over greater distances than 
the average. On the other hand we also know from ex- 
periments B and simulations that on the long time 
scales of diffusion the NGP drops, which indicates that 
the heterogeneity decays. 

Starting from these two points, we build a simple 
model to explain the time dependence of the NGP. We 
make three assumptions, all based on previously known 
results. First, all atoms have a vibrational mean square 
displacement, increasing with temperature cx T and giv- 
ing an initial ct2(t) « 0.2. Secondly, there are groups of 
mobile atoms which jump collectively. And thirdly, af- 
ter such a jump some atoms can leave while others enter 
a mobile group. Fig. 3 depicts such a collective jump 
schematically. 




FIG. 3. Schematic representation of a collective jump (left 
configuration to right configuration). Mobile and immobile 
atoms are indicated by full and empty circles, respectively. 
The grey circles show the atoms which have changed their 
group. After the jump one previously mobile atom has left 
the group, i.e. it will not participate in the next jump. An- 
other atom has joined the group instead. 

In order to make this model more tractable we use the 
most simple approximations. We have, however, checked 
carefully that the results do not depend on these de- 
tails. We consider a system formed of several groups, 
each group containing 10 atoms. Each of these groups 
of atoms can jump collectively over some barrier into an 
adjacent configurational minimum position of the under- 
lying energy landscape J|^-|35| . For simplicity we take a 
constant distribution of the activation energies. Further- 
more, we take a constant probability for jump reversal, 
ranging from 0.1 to 0.5. A backward jump means that the 
group of atoms returns to the previous positions. Each 
time a group of atoms jumps, one atom of the group 
is exchanged with an atom of another randomly chosen 
group. This accounts approximately for the changed lo- 
cal environment of the atoms after the jump which will 
cause some atoms to be in more stable neighborhoods 
and some others to become more unstable in turn. In 
the long time limit, this exchange procedure leads to the 
reduction of the heterogeneity since each atom will ulti- 
mately participate in the diffusion. To be effective, an en- 
ergy barrier has to have a minimum hight, depending on 
temperature. If the system is at a sufficiently high tem- 
perature the hopping over the lowest barriers will merge 



with the anharmonic vibrations. The atoms will instead 
be restrained by the next barrier of sufficient hight. Seen 
the other way round, when the system is cooled down, 
it becomes affected by more and more fine details of the 
energy landscape. To take this into account, we define 
at each temperature a minimum hight for a barrier to 
be effective for relaxations, -Emm = ksT. Lower barri- 
ers are replaced randomly by barriers above this thresh- 
hold, keeping a flat distribution of barrier hights. Finally 
the vibrations are accounted for by an "instantaneous" 
mean square displacement of the atoms. With this simple 
model in mind, we perform Monte Carlo-like simulations. 




t (time units) 

FIG. 4. Non-Gaussianity parameter multiplied by time ver- 
sus time, computed for a simple model of relaxation (see text) 
at different temperatures in the super-cooled liquid and in the 
glass. 

We build a system of 2000 atoms, divided into 200 
groups having different activation energies. At each time 
step we check for each group i whether this group will 
jump with the associate probability exp (—Ei/ksT). If 
the group jumps, we move all the atoms of this group, 
and then exchange one atom of the group with another 
atom. We have checked that neither the exchange rate 
nor the probability of jump reversal change the exponent 
of the power law. They merely affect the absolute value 
of the non-Gaussianity parameter. We repeat this pro- 
cedure for several time steps and starting configurations. 
During these simulations the non-Gaussianity parameter 
is computed at different temperatures from 0.42, which 
in our units corresponds to the liquid, to 0.14 which is 
under an equivalent T g . Fig. 4 shows the result of these 
simulations. 

As one can see this simple model reproduces the power 
law oc2(t) oc y/i, found in the molecular dynamics simula- 
tion, which is striking due to the simplicity of the model 
and even more so to the fact that all details of the in- 
teraction in the materials are neglected. In other words, 
the model can be applied to any kind of glass-forming 
material, once one assumes that the relaxations are by 
collective hopping of groups of particles. Therefore, we 
think that the behavior of the non-Gaussianity parame- 
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ter is universal and will also be observed in oxide glasses 
such as silica and in polymeric glasses. 

In this work we have concentrated on the intermediate 
time regime which is more easily accessible to experi- 
ment. There remain some other interesting questions. 
Comparing Fig. 2 with the predictions of the trapping 
diffusion model |56| we see a difference for long times. 
The decay of the NGP in the simulation is clearly slower 
than the predicted X/t. The reason for this is not yet 
understood. Another open question is the NGP at tem- 
peratures much below T g . Lowering the temperature the 
spectrum of the activation energies should eventually be- 
come important. 

To conclude we have presented the results of two inde- 
pendent molecular dynamics simulations on completely 
different systems and of a simple model. All these re- 
sults show the same power law for the non-Gaussianity 
in the intermediate time range, corresponding to the (3- 
relaxations in under-cooled liquids and in glasses. This 
increase of the non-Gaussianity, i.e. of the heterogeneity 
of the relaxations, proportional to y/i, can be understood 
as resulting from the collective hopping of groups of par- 
ticles. Assuming that this mechanism is common to all 
kinds of glass-formers, we believe that heterogeneity will 
always increase in the intermediate time regime domain 
following the power law y/t, at all temperatures and in 
all kinds of materials. 
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